/* 
 * File:   Raiz.cpp
 * Author: jorgeorm
 * 
 * Created on 25 de diciembre de 2011, 12:41 PM
 */

#include "Raiz.h"
#include <cstdio>
#include <iostream>
#include <cmath>
using namespace std;

Raiz::Raiz() {
    this->cotaError=0.005;
}

Raiz::Raiz(long double cotaE)
{
    this->cotaError=cotaE;
}

//Raiz::Raiz(const Raiz& orig) {
//    this->cotaErrorBiseccion=orig->getCotaErrorBiseccion();
//}

Raiz::~Raiz() {
    this->cotaError=0;
}

long double Raiz::biseccion(long double  xinf, long double xsup){
    int cotaIter =0; // Numero de iteraciones a realizar para encontrar la mejor raiz
    long double raiz=0; //Valor de la raiz actual
    long double raizAnterior=0; //Valor de la raiz calculada anteriormente
    long double limSup=xsup; //Limite superior del intervalo en el que se encuentra operando
    long double limInf=xinf; //limite inferior del intervaloi en el que se encuentra operando
    long double fRaiz=0; //Función evaluada en el valor actual de la raiz
    long double flimInf=funcion(limInf); //función evaluada en el limite inferior
    long double errorAproxPercnt=abs(limSup-limInf/limSup+limInf)*100.0; //error aproximado de la iteración actual
//    long double errorAproxAnterior=INT_MAX; //Error aproximado de la iteración anterior
    
    //Calculo de numero de iteraciones a realizar
    cotaIter = (int)ceil(log10((limSup-limInf)/(cotaError/100))/log10(2));
    cout<<"Cota en Iteraciones: "<<cotaIter<<endl;
    
    int contIter=0; //Contador de iteraciones desde la ultima mejora al valor del error aproximado.
    
    while (contIter<cotaIter&&cotaError<errorAproxPercnt){
        raizAnterior=raiz;
        raiz=(limSup+limInf)/2;
        fRaiz=funcion(raiz);
        contIter++;
        if(fRaiz!=0) {
//            errorAproxAnterior=errorAproxPercnt;
            errorAproxPercnt=(abs((raiz-raizAnterior)/raiz)*100.0);
        }
        else {
            errorAproxPercnt=0;
            break;
        }
//        Reseteo el contador en caso de haber mejora
//        if(errorAproxPercnt<errorAproxAnterior){
//            contIter=0;
//        }
        
//        Prueba de teorema de Bolzano
//        
        long double bolzano=fRaiz*flimInf;
        if(bolzano<0)limSup=raiz;
        else if (bolzano>0){
            limInf=raiz;
            flimInf=fRaiz;
        }
    }
//    Logs
//    
    printf("Raiz: %#8.10Le ±Δ= %#8.10Le\n", raiz, abs(raiz-raizAnterior));
    printf("f(Raiz)= : %#8.10Le \n", fRaiz);
//    cout<<"Raiz: "<< raiz <<" ±Δ= "<<abs(raiz-raizAnterior)<<endl;    
    cout<<"Cota superior del error: "<<this->cotaError<<"%"<<endl;
    cout<<"Error Aproximado: "<< errorAproxPercnt<<"%"<<endl;
    cout<<"Iteración No: "<< contIter <<endl;
    
    return raiz;
}

long double Raiz::newtonRaphson(long double xInicial, int cotaIter){
    long double errorAproxPercnt=100;
    long double raiz=xInicial;
    long double raizAnterior=0;
    long double fRaiz=funcion(raiz);
    long double fDeriv=derivadaFuncion(raiz);
//    long double raizReal=inversa()
    int contIter=0;
    cout<<"ContIter= "<<contIter<<" CotaIter: "<<cotaIter<<" cotaError: "<< cotaError<<" errorAprox: "<<errorAproxPercnt<<endl;
    while (contIter<cotaIter&&cotaError<errorAproxPercnt){
        raizAnterior=raiz;
        if(fDeriv!=0)
        {
            raiz=raizAnterior-(fRaiz/fDeriv);
            printf("Iter: %d, Raiz: %#8.10Le ", contIter,  raiz);
            

        }
        else {
            cout<<"Hubo división por cero, luego no se puede proyectar a un corte con el ejex"<<endl;
            break;
        }
        
        fRaiz=funcion(raiz);
        printf("FRaiz: %#8.10Le ", fRaiz);
        fDeriv=derivadaFuncion(raiz);
                
        contIter++;
        if(fRaiz!=0) {
//            errorAproxAnterior=errorAproxPercnt;
            errorAproxPercnt=(abs((raiz-raizAnterior)/raiz)*100.0);
//            errorAprox=(-1.0*(derivadaDosFuncion(raiz))/(2*fDeriv))*pow(errorAprox,2.0);
            
        }
        else {
            errorAproxPercnt=0;
            break;
        }
                printf("ErrorAprox: %#8.10Le %s\n", errorAproxPercnt,"%");

    }
    //    Logs
//    
    cout<<"----------------------------FIN----------------------------"<<endl;
    printf("Raiz: %#8.10Le\n", raiz);
    printf("f(Raiz)= : %#8.10Le \n", fRaiz);
//    cout<<"Raiz: "<< raiz <<" ±Δ= "<<abs(raiz-raizAnterior)<<endl;    
    cout<<"Cota superior del error: "<<this->cotaError<<"%"<<endl;
    cout<<"Error Aproximado: "<< errorAproxPercnt<<"%"<<endl;
    cout<<"Iteración No: "<< contIter <<endl;
    
    return raiz;
}

long double Raiz::muller(int cotaIter, long double xInicial, long double hInicial){
    long double x0, x1, x2, dr, a, b, c, h0, h1, d0, d1, fx2, fx1, fx0, errorAprox,fraiz;
    x2=xInicial;
    x1=xInicial+hInicial;
    x0=xInicial-hInicial;
    fx2=funcion(x2);
    fx1=funcion(x1);
    fx0=funcion(x0);
    errorAprox=100;
    
    long double raiz=xInicial;
    long double condiError=cotaError*raiz;
    int contIter=0;
    long double raizCuadrad=0;
    long double divisor=0;
    do
    {
        contIter++;
        h0=x1-x0;
        h1=x2-x1;
        d0= (fx1-fx0)/h0;
        d1= (fx2-fx1)/h1;
        a=(d1-d0)/(h1+h0);
        b=(a*h1+d1);
        c=fx2;
        raizCuadrad= sqrt(pow(b,2.0) - 4*a*c);
        if(abs(b+raizCuadrad)>abs(b-raizCuadrad))divisor= b+raizCuadrad;
        else divisor= b-raizCuadrad;
        dr=-2.0*c/divisor;
        raiz=x2+dr;
        fraiz=funcion(raiz);
        if(fraiz!=0) errorAprox=abs(((raiz-x2)/raiz)*100.0);
        else errorAprox=0;
        cout<<"Iteracion: "<<contIter;
        printf(" Raiz: %#8.10Le, Error Aproximado: %#8.10Le, f(Raiz)= %#8.10Le\n", raiz, errorAprox, funcion(raiz));
        
        if(errorAprox<cotaError || contIter>=cotaIter) break;
        condiError=(cotaError/100.0)*raiz;
        x0=x1;
        x1=x2;
        x2=raiz;
        fx0=fx1;
        fx1=fx2;
        fx2=funcion(x2);

//    }while(errorAprox>=cotaError && contIter<cotaIter);
    }while(abs(dr)>=condiError && contIter<cotaIter);
    cout<<"--------------------------FIN--------------------------"<<endl;
    printf("Raiz: %#8.10Le, ErrorAproximado: %#8.10Le\n", raiz, errorAprox);
    cout<<"Cota superior del error: "<<cotaError<<endl;
    printf("f(Raiz): %#8.10Le\n", funcion(raiz));
    cout<<"Iteracion: "<<contIter<<endl;
    
    return raiz;
}

long double Raiz::secante(long double  x1, long double x2, int cotaIter){
    long double errorAproxPercnt=100;
    long double valorxi=x1;
    long double valorxi_1=x2;
    long double fvalori=funcion(valorxi);
    long double fvalori_1=funcion(valorxi_1);
//    long double raizAnterior=0;
    long double divisor=0;
//    long double raizReal=inversa()
    int contIter=0;
    while (contIter<cotaIter&&cotaError<errorAproxPercnt){
        long double raiz=0;

//        raizAnterior=raiz;
        divisor=fvalori_1-fvalori;
        if(divisor!=0)
        {
            raiz=valorxi-(fvalori*(valorxi_1-valorxi)/divisor);
        }
        else {
            cout<<"Hubo división por cero, luego no se puede proyectar a un corte con el ejex"<<endl;
            break;
        }
        
        fvalori_1=fvalori;
        fvalori=funcion(raiz);
        
        valorxi_1=valorxi;
        valorxi=raiz;
                
        contIter++;
        if(fvalori!=0) {
//            errorAproxAnterior=errorAproxPercnt;
            errorAproxPercnt=(abs((valorxi-valorxi_1)/valorxi)*100.0);
//            errorAprox=(-1.0*(derivadaDosFuncion(raiz))/(2*fDeriv))*pow(errorAprox,2.0);
            
        }
        else {
            errorAproxPercnt=0;
            break;
        }
    }
    //    Logs
//    
    printf("Raiz: %#8.10Le\n", valorxi);
    printf("f(Raiz)= : %#8.10Le \n", fvalori);
//    cout<<"Raiz: "<< raiz <<" ±Δ= "<<abs(raiz-raizAnterior)<<endl;    
    cout<<"Cota superior del error: "<<this->cotaError<<"%"<<endl;
    cout<<"Error Aproximado: "<< errorAproxPercnt<<"%"<<endl;
    cout<<"Iteración No: "<< contIter <<endl;
    
    return valorxi;
}

//Seters and Getters 
//
void Raiz::setCotaError(long double cota){
    this->cotaError=cota;
}

long double Raiz::getCotaError(){
    return this->cotaError;
}


